Modification of the Gay-Berne potential for improved accuracy and speed 
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A modification of the Gay-Berne potential is proposed which is about 10% to 20% 
more speed efficient (that is, the original potential runs 15% to 25% slower, depending 
on architecture) and statistically more accurate in reproducing the energy of inter- 
action of two linear Lennard- Jones tetratomics when averaged over all orientations. 
For the special cases of end-to-end and side-by-side configurations, the new potential 
is equivalent to the Gay-Berne one. 



'Electronic mail: 



rasmus.persson@chem.gu.se 



The Gay-Berne potential is a coarse-grained potential designed specifically to reproduce 
the energy of interaction between two linear, tetratomic Lennard- Jones molecules. As such, 
it allows the approximate treatment of a linear, rigid polyatomic segment by use of a single 
interaction site, reducing the amount of looping - and hence the calculational effort - drasti- 
cally when calculating energies, forces and torques in molecular mechanics. The Gay-Berne 
potential is not the only single-site potential for linear molecules or molecular moieties lack- 
ing substantial electric multipole momomentsj^ but it is arguably the most widely used. 
In this Note, we present a variation of the Gay-Berne potential which undeniably improves 
even further on the computational economy and also on the accuracy. For ease of reference 
and of comparison, the original functions are also reproduced. 

The Gay-Berne potential, for intermolecular distance r, is given in eq. ([I]). 
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In this equation the e(u\,u 2 ,r) and a(ui, u 2 , r) are auxiliary functions that depend both on 
the molecular orientations, expressed through the unit vectors Si and u 2 , and on the orien- 
tation of the intermolecular vector, expressed by r . Gay and Berne^ provide the following 
forms for these functions: 
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where do and eo are constants, and \ and x' relate to the anisotropy of the interaction and 
are defined by x = (J 2 ~ l)/(^ 2 + 1) an d x' = (v^ — l)/(\/rf + 1), where / is the steric 
length-to-breadth ratio (the axial ratio) and d is the ratio of the binding energy in the side- 
by-side configuration to that of the end-to-end configuration. A comparison of the energy 
of interaction between the model system it was designed to replicate, and the Gay-Berne 
potential itself is shown in Figure [TJ using the original parameter set advocated by Gay 
and Berne themselves. We shall keep the general form of eq. ([T]), but replace the auxiliary 
functions a and e by more efficient alternatives. 

The choices for a and e that we make are inspired by Ref. [TUl in which the equipotential 
surface of the Lennard- Jones potential was redefined to that of an ellipsoid. There, the 
philosophy of Berne and Pechukas 6 and of Kihara^ - that the anisotropic intermolecular 
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FIG. 1. Energy of interaction of two linear, tetratomic Lennard- Jones molecules of bond length 
2cr/3 (full line) and their Gay-Berne approximation (dashed line) with I = 3 and d = 5. The 
binding energy of the side-by-side configuration is set to unity in both cases. The unit of length is 
a of the Lennard-Jones interaction. 

potential is obtained by replacing the distance scale by an anisotropic function - was strictly 
adhered to. In abandoning that principle for the prescription of Gay and BerneP so that 
we displace rather than dilate the interaction potential, using the anisotropic a function of 
Ref. [Eland its close analog for the e function in eq ([T]), we arrive at the modified Gay-Berne 
potential. To avoid a notational cluttering and as confusion is unlikely to arise, we do not 
distinguish these new functions through notation from their Gay-Berne counterparts. The 
new angular functions are represented by 
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These new functions give the same dashed curves in Figure [T] as the original ones. Important 
as such agreement might be, those two geometries only represent a subset of all possibil- 
ities; a more thorough test of the model is to calculate its average error over all possible 
configurations. The average absolute error, denned by ( Inexact /-^modei — 1|) where the angle 
brackets denote an unweighted angular average, is shown in Table [I] for both the original 



and the modified Gay-Berne potential as a function of intermolecular separation. E exact is 
the energy of the dimer of the linear tetratomic, as used in the original paper^ and -E mo dci 
is the energy of interaction of the Gay-Berne potential, original or modified, in the same 
geometry. Both energies are normalized by their respective side-by-side binding energy. It 
is worth pointing out that this comparison is made with the original Gay-Berne parameter 
set, optimized for this situation, and that further improvement is perhaps possible with ad- 
justments. It is a chief strength of the modification that it does not need reparametrization 
to improve upon the original. 

TABLE I. Angular- averaged unsigned errors for selected intermolecular distances r for the original 
and modified Gay-Berne potential with respect to the linear Lennard-Jones tetratomic, obtained 
by unbiased Monte Carlo sampling over 10 8 cycles. Identical values of I = 3 and d = 5 are used in 
the comparison. 
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One further virtue of these new angular functions is the ease with which they can be 
extended to cover also the case of the interaction between non-identical elongations and 
binding energies. The natural generalization of the anisotropic functions is 
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where subscripts 1 and 2 distinguish between two different types of bodies. For a complete 
specification of the unlike interaction, also the corresponding eo and <j constants need be 
specified, something which can be done by the standard Lorentz-Berthelot combining rules. 
Corresponding generalizations for the original Gay-Berne potential have been presented by 
Cleaver et alW and by Berardi et a/P^ The latter of these generalizations covers also the 
case of biaxial molecules or molecular segments. 

In closing, we make a quick comparison of the computational complexity of the modified 
and original potential. It is clear that eqs ^ and ^ require a much greater number 
of arithmetic operations than do eqs Q and ([5]). Eqs ^ and ^ also require the non- 
elementary square-root operation. In contrast, eqs Q and ^ require only the elementary 
arithmetic operations multiplication and addition. In the author's hands, the subroutine for 
the original potential requires approximately 15% to 25% more computer processing time 
(tested on two different computers). Using fast floating-point (e. g. less accurate) libraries 
is likely to diminish this difference somewhat. In actual simulations, the speed benefit will 
be somewhat less depending on system size. 
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